First verification of generic fidelity recovery in a dynamical system 
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We study the time evolution of fidelity in a dynamical many body system, namely a kicked Ising 
model, modified to allow for a time reversal invariance breaking. We find good agreement with the 
random matrix predictions in the realm of strong perturbations. In particular for the time-reversal 
symmetry breaking case the predicted revival at Heisenberg time is clearly seen. 
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I. INTRODUCTION 

The random matrix (RMT) description of fidelity de- 
cay pj has been successful in describing the behavior of 
dynamical models and actual experiments 0, El • While 
in this early work linear response (LR) results were ex- 
ponentiated (ELR) to describe long-time behavior quite 
well, in a more recent paper it was shown, that the RMT 
model presented in p| can be solved exactly in the limit 
of large Hilbert space dimension 0, |j| . These exact re- 
sults show a puzzling, though in absolute terms weak, 
revival of the ensemble averaged fidelity amplitude at the 
Heisenberg time. The predicted revival is strongest for 
the Gaussian symplectic ensemble (GSE), and weakest 
for the Gaussian orthogonal ensemble (GOE). For the 
time being experiments are far from the realm of param- 
eters where this effect or any significant deviations from 
ELR can be seen; neither have such deviations been seen 
in dynamical models. 

In the present paper we intend to show that this ef- 
fect exists in a dynamical many particle model in which 
we can carry calculations beyond the Heisenberg time for 
very large systems, and for large perturbation strength. 
The system we shall consider is a modification of the 
kicked Ising chain [|| . Since we are dealing with a kicked 
system we have to focus on the Floquet operator and 
hence we shall have to compare to the so called cir- 
cular ensembles of unitary matrices. Note that in the 
large dimension limit the statistics for spectral fluctua- 
tions coincide. We introduce a multiple kick in order to 
break an anti-unitary symmetry that plays the role of 
time inversion and which we still call time reversal in- 
variance (TRI). To consider TRI breaking evolution is 
useful as the revival at Heisenberg time is in absolute 
terms three orders of magnitude larger for the Gaussian 
unitary ensemble (GUE) than for the GOE. The unper- 
turbed Hamiltonians used for both the TRI and the non- 
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TRI situations are fixed such as to yield very good agree- 
ment with the spectral statistics expected for the GOE 
or GUE, respectively. This choice is certainly not unique 
as we could further modify the kick-structure, but this 
is not contemplated. The perturbation strength, that 
enters the RMT model as a free parameter, can be deter- 
mined from the dynamics via correlation functions. In 
both cases we shall thus find good agreement between 
the fidelity of the dynamical model and RMT without 
any free parameters. 

The effects we look at are of the order of a part in a 
thousand or less for the ensemble average, and we wish 
to analyze whether they can be seen in an individual 
system. For this to be the case we need large enough 
Hilbert spaces such that the time evolution samples a 
large domain of eigen-frequencies of the system, which 
in turn gives the opportunity for self-averaging. In the 
specific case of kicked spin chains we can numerically 
handle time evolution for up to 24 qubits for moderate 
times. 

Fidelity describes the evolution of a cross-correlation 
function of a given initial state under two differ- 
ent Hamiltonians or, equivalently, the evolution of 
the auto-correlation function under the so-called echo- 
dynamics 0- The meaning of the latter term is a com- 
position of the forward time evolution by some Hamil- 
tonian and backward time evolution by a slightly differ- 
ent one. For an initial state \tp(0)) evolving under some 
unperturbed unitary dynamics Uo(t) and the same state 
evolving under a perturbed but still unitary dynamics 
U e (t) the fidelity amplitude reads as 

f(t) = {i>(o)\uhWe(t)\m)- (i) 

Fidelity itself is defined as F(t) = \f(t)\ 2 . 

We shall first review the RMT results, and give a 
brief resume of the linear response and exact solutions. 
Next we discuss the generalized kicked spin model and 
show that we can reach sufficiently large Hilbert spaces 
and sufficiently strong perturbations to be in the desired 
regime. First we check that for a suitable non-integrable 
parameter regime of the model we can have GOE statis- 
tics for the TRI conserving model and GUE statistics for 
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the TRI violating model. Finally, we compare the model 
results for fidelity decay with the corresponding RMT 
predictions. 



II. RANDOM MATRIX RESULTS 

Let Hq be the unperturbed Hamiltonian, chosen from 
one of the Gaussian ensembles, and 



(2) 



the perturbed one. This somewhat unusual definition 
of the perturbation strength e has been applied for con- 
venience of comparison with the RMT results. The fi- 
delity amplitude is then given by Eq.Q with U{t) = 
cxp(—2iriHot) and U e (t) = cxp(—2iriH e t). It is assumed 
that Hq has mean level spacing of one, and thus t is given 
in units of the Heisenberg time th ■ The variance of the 
off-diagonal elements of V is chosen to be one. 

Gorin et al calculated the Gaussian average of the 
fidelity amplitude in the regime of small perturbations 
using the linear-response approximation, 



/(f)~l-eC(t). 



where C{i) is given by 
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1 — b2,p(t) is the spectral form factor, and (3 is the uni- 
versality index, i. e. (3 = 1 for GOE, /3 = 2 for GUE, and 
/3 = 4 for GSE. For an explicit calculation, knowledge of 
the spectral form factor is thus sufficient. Using the ELR 
approximation, 
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= c(t) 
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the authors were able to describe quantitatively the cross- 
over from Gaussian to exponential decay with increasing 
perturbation strength. 

It is obvious that the linear-response approximation 
must break down for large perturbations. In a recent 
paper Stockmann and Schafer applied super-symmetry 
techniques to calculate Gaussian averages of the fidelity 
decay for arbitrary perturbation strengths In or- 
der to apply the super-symmetry framework, the authors 
considered the state averaged version of Eq. (JIJ using the 
Hamiltonians define above and the perturbation param- 
eter e 
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e 2mH c t -2mH t\ 
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where N denotes Hilbert space dimension, hence making 
it independent of the state considered. 

The exact results obtained in Ref. |f| are in good agree- 
ment with the linear-response result for very small per- 
turbations and also confirm the validity of the ELR ap- 
proximation for moderate perturbation strengths, and 
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FIG. 1: We show the probability density for the nearest 
neighbor spacing for 18 qubits and averaged over the dif- 
ferent relevant TLk spaces (k — 1,...,8). Numerical results 
are shown in green circles for the model corresponding to the 
GUE and red triangles for the GOE. Parameters are M = 2, 
b {1) = (1, 1,0), and b (2) = (1.4,0, 1.4) for the GUE case and 
= (1.4, 1.4, 0), for the GOE case. We plot the exact RMT 
result for both the GOE (dashed red line) and the GUE (con- 
tinuous green line) since comparison with Wigner surmise is 
not satisfactory. 



indeed whenever the fidelity is not too small. However, 
these calculations also revealed an important generic fea- 
ture, which motivated the present work: for very strong 
perturbations, the fidelity amplitude shows a partial re- 
covery at the Heisenberg time. 

For the GUE case the result for the fidelity amplitude 
reads 



/(*) 
where 
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s(x) 



sinh(a;) 



(7) 
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and s'(x) denotes its derivative. For the GOE case, the 
result is not so simple and is expressed in terms of one- 
dimensional integrals [5|- 



III. THE MULTIPLY KICKED ISING MODEL 

We next introduce the multiply kicked Ising model 
(MKI) as a modification of the kicked Ising model in- 
troduced in 6] . The evolution operator corresponding to 
one time step is 



M 

UmKI = J\ U lsing U k f ck 



(9) 
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FIG. 2: Same parameters and sign conventions as in Fig. Q 
We show the spectral form factor averaged over all relevant 
Hk's. In order to smooth the function, we also average over 
a window in t of size 0.154h- 



with 



L-l 



exp -*J$>K+i > 



3=0 
L-l 



^kick 



exp (-z^ &(")•«?,], 

3=0 



(10) 
(11) 

and periodic boundary conditions ctl = &q- Here a^' v,z 
are the Pauli matrices corresponding to spin 1/2 at po- 
sition j, and <7j = (ctJ, crj, cr|). The system defined in 
Eq.® represents a periodic 1-d array of L spin 1/2 par- 
ticles which periodically receive a sequence of M differ- 
ent kicks of instantaneous magnetic field pulses equally 
spaced in time. The free evolution (f7i s ; ng ) is simply given 
by the Ising interaction between nearest neighbors, with 
dimensionless strength J which we fix to 1 in the numer- 
ical simulations. Each of the M kicks is an instantaneous 
and homogeneous magnetic field characterized by the di- 
mensionless vector . Due to the two body interaction 
of the system, we are able to evaluate efficiently the evo- 
lution of any state 0, and hence, as noted in Q, also its 
spectrum. In the previous model 6] the same kick was 
applied periodically, that is M = 1. 

Since we wish to compare the behavior of fidelity in 
our model-system with the average behavior of an en- 
semble of random matrices, it is necessary to identify the 
correct symmetry class. To this end, we have to study 
the symmetries of the MKI model. As couplings between 
neighboring qubits as well as the effect of the kick on all 
qubits are equal, it is easy to see that we formally have a 
ring of qubits and thus a rotational as well as a reflection 
symmetry. 

Let us define the rotational operator R over the ele- 
ments of the computational basis, as R\ioii . . . = 
|*l-i*o • • ■ *L-2), where each ij € {0, 1}, and extend it to 
the other elements of the Hilbert space, requiring linear- 



ity. If we visualize our model as a ring, the action of this 
operator is to rotate the particles in the ring by one posi- 
tion. This operator defines a rotation group Cl- We can 
further define the reflection operator P on the compu- 
tational basis as P\ioi\ . . . il-i) = |*l-i*l-2 • ■ ■ *o) & n d 
again extend it to the whole space requiring linearity. 
The two operators do not commute and form the group 
Cl,v of rotations and reflections of a ring of L elements. 

Note that the [Umki,P] = [Umki,R] = 0, hence the 
eigenvalues, exp(27rj&/L) with k E Z/i, of the rotation 
operator define the invariant subspaces which are degen- 
erate with two parities except for the cases k = and 
k = L/2 the latter for even L only. In the latter cases 
invariant subspaces of well defined parity can be defined 
depending on the behavior under reflections, though for 
a given set of qubits not necessarily both parities need 
occur. The dimension for any set with fixed k is approx- 
imately equal to 2 L /L. 

Since we do not identify other symmetries in the sys- 
tem, we expect to find a range of parameters in which 
each of the blocks of C/mki corresponding to the min- 
imal invariant subspaces, behave as a typical member 
of the GUE. This was numerically checked by analyzing 
the spectra in each of the invariant subspaces, for the 
parameters shown in the caption of Fig. ^ The result- 
ing nearest neighbor statistics and spectral form factors 
(both defined in 0) are presented in Fig. ^ and Fig. 
respectively. Other statistical tests like the number vari- 
ance, skewness and excess were also applied with excel- 
lent agreement with the expected behavior (not shown). 
The success of these tests is a strong indicator that in- 
deed no other symmetries are present. 

Consider now the case in which we have only one pe- 
riodic kick, that is M = 1 in E/mki- Rotating each in- 
dividual spin, the magnetic field (b^) can be made to 
have only components in the xz plane. The Hamilto- 
nian will have only real components in the basis in which 
<tJ' z are real, hence an anti-unitary symmetry becomes 
evident. This symmetry K, (complex conjugation in the 
basis mentioned above) is not time reversal, since spin is 
being reflected in the xz plane of each qubit, instead of 
being reversed. Furthermore, this symmetry will change 
the sign of k and hence, by itself, will not provide an anti- 
unitary symmetry within the invariant subspaces. Re- 
calling that P also reverses the sign of k, we observe that 
the combined operator P)C, will provide each of the in- 
variant subspaces with an anti-unitary symmetry, which 
we call time reversal invariance (TRI). We expect then 
that for M = 1 the system will have a range of parame- 
ters in which each invariant subspace will behave as a typ- 
ical member of the GOE for some values of the parame- 
ters which are sufficiently well separated from the exactl 
solvable cases of longitudinal and transverse fields [l 
Numerical evidence favoring this statement is presented 
in Figs. [U and 

Effectively we will have, for each of the invariant sub- 
spaces a Heisenberg time given by in ~ 2 L / L (the ap- 
proximate number of levels) except for TRI subspaces for 
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which we shall have in ~ 2 L ^ 1 /L. 

The perturbation we shall consider is a kick in the x 
direction: 

A=$>* (12) 

j=o 

in the Floquet propagator [Eq.©] 

Umki,s = t/MKiexp(-2<L4). (13) 

Comparing the linear response formula for the dynamical 
model 0, 

/(*) = !- If E (14) 
t'=-t+i 

with the correlation function C(i) = 2-- L TrA[/ t (t)AZ7(i) 
with the random matrix model, we connect the pertur- 
bation strengths 

e = 2 L 6 2 a . (15) 

Here a = lim^oo lim^^oo i ^*, = t C(t') is the inte- 
grated correlation function of the unperturbed dynamics. 

Note that the perturbation will have all the symme- 
tries presented, so that it will not mix the different in- 
variant subspaces. Thus, if we consider an initial state 
with components from all the Tik spaces, we will obtain 
an average of fidelity over L different initial conditions 
each in a different 7ik space. 



IV. RESULTS 

Since the effect we want to observe is extremely small, 
we must have a full understanding of the most impor- 
tant causes of deviation in the evaluation of fidelity in 
Eq. - Note that here we are taking two averages, one 
of them over the ensemble of Hamiltonians and the other 
one over the ensemble of initial conditions. The first of 
this averages we are not going to evaluate in our nu- 
merical model since we want to show that an individual 
quantum chaotic system (in the limit of large Hilbcrt 
space dimension) actually behaves according to the ran- 
dom matrix ensemble average over the appropriate sym- 
metry class. The other average (over the Hilbert space) 
is possible to evaluate in an exact way, taking the average 
over an orthonormal basis resulting in a trace operation. 
In practice, for very large Hilbert spaces this method 
is not efficient. Instead we shall use an approximation 
which introduces an error that can be made as small as 
desired. 

Let us first comment on the latter source of error. For 
any operator A, we have that 

TrA = (A), , (16) 
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FIG. 3: Different perturbations for the GUE system. We 
have here 16 qubits and 10 initial conditions, in Fig. |S] Dots 
show the numerical calculation (real part) and thin curves 
show RMT prediction whereas red dashed curves show the 
ELR approximation. We set perturbations to ei = 5.15 (in 
this case all the three curves are overlapping), 62 = 10.3, 
e 3 = 15.455, e 4 = 20.6. 

where (•) . denotes the average over Gaussian random 

states: = X^i^iK); where {\i)}i=i,...,N is an or- 
thonormal basis and independent complex Gaus- 
sian distributed random numbers with standard devia- 
tion 1 / >/~N. Then, we have that 

Of course letting m — > oo will make the average exact. 
But for large N it turns out that already taking a small 
number of states is enough to have a very good approxi- 
mation of the expected value of A. So we shall evaluate 
Eq. iJBJ with the aid of Eq. l(T7l) . In particular for the 
sizes of the Hilbert spaces considered here (upto 2 20 ) the 
number of states needed to achieve the precision required 
is much smaller than the dimension of the space. In our 
case we want to evaluate the expectation value of the 
echo operator 

M t = u\{t)U{t) (18) 

for different values of t. For a fixed value of t we have a 
distribution of values of (t/j\M t \ijj) and hence an associ- 
ated standard deviation. Although at t = this standard 
deviation is zero, after a transient time it approaches an 
stationary value, determined by finite size effects. We 
shall call this value Cfmite average- Then the error in the 
evaluation of TrM t will be Cfinitc average/ y/rn if we take m 
sample states. Estimating numerically this value is triv- 
ial since for m different realizations we can evaluate the 
standard deviation at each time and then average over 
time, for t's larger than the transient time. 

Considering only one particular Hamiltonian will also 
cause some deviations from the exact RMT formula: re- 
call that to obtain, say Eq.Q), we averaged over an en- 
semble of Hamiltonians. Even if we obtain the exact trace 
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FIG. 4: Different number of qubits for a same effective per- 
turbation for the GUE system. The perturbation strength is 
e = 20. We vary the number of qubits, red triangles being 
10, green squares 12, and blue circles 16. We have, in all 
cases, 200 initial conditions. We only show the points where 
Re/ > 0. The red dashed line is the ELR prediction. 



of the echo operator its value will fluctuate in time around 
the ensemble averaged value. These fluctuations can be 
characterized with a standard deviation <7intrinsic- For a 
fixed Hamiltonian the only way to decrease this number 
is to increase the dimension of the system. In practice, 
for computational reasons, we can only achieve results 
for up to 20 qubits for times of the order of Heisenberg 
time. 

Since those two effects can be assumed to be com- 
pletely independent, the total average deviation of Eq.© 
will be given by 



""total "intrinsic 



""finite average 



(19) 



We can estimate a posteriori this quantity very easily, 
since we know that the imaginary part of the ensemble 
averaged fidelity amplitude is zero for all times: Im/(i) = 
0. Hence computing the fluctuations of the imaginary 
part of the data obtained will give us ertotai and hence 
^intrinsic (recall that evaluating a^ nltc average is trivial). 
The knowledge of extrinsic will give us an estimate for the 
minimal value of m needed to estimate the value of TrM t 
with sufficient accuracy. In the results shown in this sec- 
tion, we always have <7 nni tc average l\frn < ^intrinsic. 

We now turn our attention to the numerical results 
obtained. We shall set &W = (1.4, 1.4,0). For the GOE 
case (M = 1) wc do not need to specify any other param- 
eter except for the number of qubits and the strength of 
the perturbation. The resulting value of the integrated 
correlation function is found to be a ~ 1.27. In case we 
want to observe the GUE prediction we set M = 2, take 
the same as in GOE case and b {2) = (1,0,1). The 
integrated correlation function is found to be a 0.93. 

In Fig.[3]we show, for the GUE- type system, how vary- 
ing the perturbation parameter changes the behavior of 



FIG. 5: Revival for 20 qubits. The perturbation strength is 
e = 31.78. The red triangles indicate the real part, whereas 
the blue circles indicate the imaginary part (whenever they 
are greater than 0), in time steps of 525, calculated with 15 
initial conditions. The black thin curve is the RMT theo- 
retical prediction, the red dashed line is the ELR prediction 
and the gray curves indicate the RMT curve plus/minus the 
calculated intrinsic deviation. 



fidelity. Even for fairly large perturbations (ei = 5.15) 
we obtain a fairly good approximation with ELR, which 
almost coincides with the exact RMT result. For larger 
perturbations, deviation of the ELR from exact RMT are 
big enough to be observed, and indeed they are observed 
clearly for ei = 10.3. For £4 = 20.6, where we should 
observe a recovery, the behavior at those low fidelities 
is shadowed by the intrinsic error. Here we chose suffi- 
ciently many initial conditions, so that the finite average 
error is considerably smaller than the intrinsic error. 

In Fig. 0| we show, again for the GUE-type system, 
how increasing the number of qubits decreases the devi- 
ation from RMT prediction. We fix the effective pertur- 
bation [which scales with the number of qubits as shown 
in Ea. H15(l ] and observe the behavior of the real part of fi- 
delity for 10, 12 and 16 qubits. Particularly for 10 qubits 
we see a noticeable deviation, probably due to the non- 
vanishing correlation function for large times. This value 
becomes smaller as the number of qubits increases but 
also we average over a larger number of invariant sub- 
spaces. In all cases we can see that the correspondence 
with exact RMT is much better than with ELR. 

In Fig. [5] we show the evolution of fidelity, for e = 31.78 
where the revival is strongest. The error due to finite 
averaging is much smaller than the intrinsic error. Wc 
cannot diminish intrinsic error further, since for 20 qubits 
we need near more than two weeks of computer time per 
initial condition, and in order to increase the system size 
by one qubit we need 4 times more CPU time for each 
initial condition, since for each time step we roughly need 
to double the number of operations and we also need to 
double the number of time steps (due to the increase 
in the Heisenberg time). We expect that for sufficiently 
large Hilbert spaces the intrinsic error will be so small 
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FIG. 6: We show that our GOE system indeed behaves as 
given by the RMT theory. Same color coding as in Fig. 
The perturbation strength is e = 10, and the number of initial 
conditions is 15. The size of the system is 20 qubits, h± — 
= 1.4. We show only the points where / > 0. The red 
dashed line is the ELR prediction. 

and the self averaging so strong that even for one initial 
condition one can clearly observe the fidelity revival. 

Finally, for completeness, we also show the behavior 
for the GOE-type system in Fig. The deviation with 
ELR is evident, and the error is within the bounds. 

V. CONCLUSIONS 

The recent exact solutions of the random matrix model 
for fidelity decay display quite unexpected characteris- 
tics for large perturbations and large times. The most 
remarkable one is the maximum that develops at Heisen- 
berg time. We showed that this phenomenon, established 



for an ensemble average, can actually be seen in indi- 
vidual dynamical systems if the Hilbert space is large 
enough: the effect is not overshadowed by noise and fi- 
nite size effects. Self-averaging is effective, and indeed we 
find good agreement with RMT results for fidelity decay 
for two individual kicked spin chain models, one with, one 
without a pseudo time reversal invariance. The results 
are achieved without any fit parameter as the pertur- 
bation strength has been determined directly from the 
model we use. 

The model used was a kicked spin chain which can re- 
produce both TRI conserving and TRI breaking dynam- 
ics. It has no obvious classical analogue, so we could feel 
reasonably confident that it would display generic behav- 
ior in a regime far from the perturbative one. Whether 
the maximum at Heisenberg time can also be expected 
for systems with a classical analog is an open and in- 
teresting question. Clearly, for strong perturbations and 
short times they should not follow RMT behavior, but 
that does not exclude, that this behavior is recovered at 
least qualitatively for times as long as the Heisenberg 
time. 
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